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ABSTRACT 

The problem of the stability of magnetic fields in stars has a long history and has been investigated in detail 
in perturbation theory. Here we consider the nonlinear evolution of a nonrotating neutron star with a purely 
poloidal magnetic field, in general relativity. We find that an instability develops in the region of the closed 
magnetic field lines and over an Alfven timescale, as predicted by perturbation theory. After the initial unstable 
growth, our evolutions show that a toroidal magnetic field component is generated, which increases until it 
is locally comparable in strength with the poloidal one. On longer timescales the system relaxes to a new 
non-axisymmetric configuration with a reorganization of the stellar structure and large- amplitude oscillations, 
mostly in the fundamental mode. We discuss the energies involved in the instability and the impact they may 
have on the phenomenology of magnetar flares and on their detectability through gravitational-wave emission. 
Subject headings: stars: neutron — gravitational waves — magnetohydrodynamics (MHD) — methods: nu- 
merical 



1. INTRODUCTION 

During at least two points within a neutron star's (NS) life, 
large-scale magnetic field rearrangement may occur. These 
are shortl y after the formation of NSs in supemovae (|Bo-| 
[nanno et a l. 2003 ), and also during the giant flares of mag- 
netars (^T hompson & Duncan,, 1996^ .Geppert & Rheinhardt] 
|2006 ). Whilst similar rearrangements may also occur in other 
stars, they are likely to be particularly significant for the 
physics of NSs, where the fields are exceptionally strong: up 
to 10^^ G at the surface of normal pulsars and 10^^ G for mag- 
netars. There are a variety of instabilities in NSs, and in a 
proto-NS mag netic fields may actua l ly have a stabil i sing e f- 
fect (see, e.g. [Miralles et al.| ( |20Q2| ); [Bonanno et aL] ( |20Q3| )), 
but we are concerned here with the fast-acting "Tayler insta- 
bility" which affects purely poloidal (or purely toroidal) mag- 
netic fields in stars. 

The magnetic-field geometry of a NS is important for the 
star's evolution, p rovides a distortion that may lead to gravi- 
tational radiation ( [Bonazzola & Gourgoulhon||1996| ), as well 
as powering the mechanisms by which these stars may be ob- 
served: the pulsar emission for normal NSs, and the X/7-ray 
emission of magnetars. It is important therefore to determine 
which models of magnetised NSs are stable equilibria. 

The study of magnetised ste llar equilibria dates back 



to 



Chandrasekhar & Fermi[ ( [1953| ). Since then, many possible 



magnetic equilibria have been studied, using both analytic and 
numerical techniques. The se have includ ed configurations 
with purely poloidal fields (|Ferraro[| 19541 Monaghan 1965; 
Bocquet et al. 199 5]) and purely toroidal fields^ (Roxburgh 
1963, ,Kiuchi & Yo'shida"2008 ), as weU as mixed poloidal- 



2008^ 



toroidal configurations (Roxburgh 1966; Has kell et al. 
Tomimura & Erig uchi[|2005f |Lander & Jones[P009[ 
etal. 2009,2010). 

However, constructing a configuration in equilibrium is 
only half the problem when modelling stellar magnetic fields; 
one also needs them to be stable over many dynamical 
timescales, since stellar magnetic fields have been observed 
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to be long-lived. This has proved to be a challenging prob- 
lem for analytic methods, which can only study the initial 
localised instability and not the resultant field configuration. 
With purely poloidal and purely toroidal fields known to be 
unstable (Markey & Tayl erll973[ [Wright|1973t |Tayler"1973) 
Flowers & Ruderman 1977]), only mixed-field configurations 
are likely to exist in stars. 

More recently it has become feasible to use numerical evo- 
lutions to study these hydromagnetic instabilities, with the 
benefit that the global behaviour of the instability may be 
studied (analytic works rely on local analyses), as well as 
the final outcome of the i nstability when th e field undergoes 
significant rearrange ment (ILander & Jones 1201 1 ; Bra ithwaite[ 
2m7|[Geppert & Rh einhardt 2006; Kiuchi et al. 201lJ). De- 
spite this recent progress, there are still very few models of 
stellar magnetic-field configurations whose stability has been 
assessed. 

The instability-induced redistribution of magnetic flux is 
potentially a very violent event and it has been suggested as a 
trigger mechanism fo r the giant flares of magnetars (Thomp-] 
son & Duncan[[199^ . This redistribution is likely to be ac- 
companied by a significant change to the mass quadrupole 
moment of a NS, making it a potentially detectable source of 
gravitational waves (GWs) ( Kashiyama & loka 201 1 Corsi & 
Owen 2011 ). For this reason it is important to understand the 
frequency, amplitude and duration these GWs may have. 

This paper is organised as follows. In Section 2 we give 
a description of our computational infrastructure and initial 
stellar models. In Section 3 we present results from our evo- 
lutions, showing the generation of an instability and the sub- 
sequent reorganisation of the magnetic field into a more stable 
configuration. We also study the GW emission from the insta- 
bility and assess its detectability. Conclusions are presented 
in Section 4. 

2. PHYSICAL SYSTEM AND NUMERICAL SETUP 

We model our initial NS as a nonrotating isolated fluid body 
in ideal magnetohydrodynamics (MHD) and with a purely 
poloidal magnetic field permeating it and extending to the ex- 
terior. The initial axisymmetric equilibrium configuration is 
generated by the LORENE code, which produces a fully rela- 
tivistic solution and consistently accounts for the metric dis- 
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Figure 1. Snapshots of the development of the instability in our fiducial star, showing projections on the (x, z) plane (upper row) and (x, y) plane (lower row) 
of the simulation at times t = 1, 3, 10 ms (left to right), respectively. Shown with vector lines are the (global) magnetic-field lines, while the colours show the 
intensity of the toroidal magnetic field only; also reported are the iso-density contours of the rest-mass density near the stellar surface. 



tortions a nd structure deformati ons induced by the magnetic 
field (see Bocquet et al. (1995)). The star is modelled as a 
poly trope with equation of state p = Kp^ , where F = 2 and 
K = 98.5, so that a NS with mass of 1.41 has a radius 
of about 12.1 km. Because the timescale for the instability is 
shorter (and computationally feasible) for stronger magnetic 
fields, we have considered stars with Bq/1{)^^ G G [2.5, 10], 
where is the value at the magnetic pole, selecting a value 
of Bq = 6.5 X 10^^ G as the reference one. 

Using these initial configurations, we perform general- 
relativistic MHD simulations in three spatial dimensions un- 
der the Cowling approximation, i.e. we do not evolve the Ein- 
stein equations but consistently solve for the MHD equations 
in a fixed and curved spacetime. This choice is motivated by 
wanting to reduce computational costs and by the fact that 
the changes in the spacetime are expected to be intrinsically 
small. As a consequence, the GW emission is computed us- 
ing the Newt onian quadrupole formula (Nagar et al. 2007; 
Baiotti et al.| [2009). The evolutions are performed with the 
WhiskyMHD code, whose properties have be en tested and 
discussed in a number of earlier papers ([Giaco mazzo & Rez- 
zolla'2007' Pollney et al.|2007[ |Giacomazzo et al. 2009||Giaj 
comazzo et al. 20 1 1 1) . Our standard n umerical setup consist s 
of a gridTwith three refinement levels dSchnetter et al.|20 Q4), 
the highest one having a resolution h /Mq = 0.17 :^ 250 m 
and covering all of the star. The outer boundaries are placed 
at a distance of 54 Mq c:^ 79 km. 

The most salient difference of the code with respect to 
the references above is in the treatment of the atmosphere. 
We recall that as customary in relativistic hydrodynamics us- 
ing finite- volume methods, we surround the star with a low- 
density "atmosphere", whose dynamics is prescribed by suit- 
able boundary conditions. More specifically, the rest-mass 
density there is set to a constant value, while the fluid ve- 



locity is reset to zero (see [Baiotti et al.| d2005| ) for details). 
Although this prescription works very well in hydrodynamic 
simulations, it becomes problematic in ideal MHD, since it 
prevents any evolution of the magnetic field in this region. For 
sufficiently strong magnetic fields this approach can rapidly 
lead to errors at the stellar surface, which prematurely ter- 
minate the simulations. To improve on our treatment of 
the atmosphere and allow for a dynamics of the magnetic 
field at surface and outside the star, we add a magnetic dif- 
fusivity term to the induction equation, which we write as 
dt{B^) = dj{v^B^ — v^B ^) + rjdj d^B^, where r] is the scalar 
resistivity (s ee [Giacomaz zo & Rez zolla (2007 ); Giacomazzo] 
|et al.| ( |2011| ) for details on the implementation in ideal-MHD 
case). Because we want to retain the ideal-MHD behaviour 
in the stellar interior and allow for an evolution in the atmo- 
sphere, we set the resistivity to zero within the bulk of the star, 
letting it increase continuously to its atmospheric value start- 
ing from a low-density region near the stellar surface. More 
specifically, we set 7/(p) = r]of{p), where p is the rest-mass 
density, f{p) is the Fermi function and r]o is a free parame- 
ter. We used a reference value of ijq/Mq =0.12, but we will 
discuss how results change for a lower value. 

Again to reduce computational costs and because we are in- 
terested in the development of magnetic-field instabilities, we 
add a perturbation designed to trigger them. From the anal- 
ysis of [Markey & Tayler] ( 1973 ) we expect a 6>-component of 
the velocity near the neutral line to induce the fastest-growing 
instability. Hence, we introduce a perturbation of azimuthal 
index m = 2 in this quantity, such that the relative change in 
the magnetic field is 10~^. We have verified that the results do 
not change qualitatively for other choices of the perturbation 
and that an instability develops even without a perturbation. 

3. RESULTS 
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Before discussing the nonlinear development of the insta- 
bility, it is useful to recall th e predictions of the perturbative 
studies of iMarkey & Tayler| (T^3 ) and |Wright| ( |1973| ) about 
its early growth. In particular, we expect that: (//the instabil- 
ity should be localised in the closed-field line region around 
the neutral line (i.e. where the poloidal magnetic field van- 
ishes); (ii) the instability should occur after about an Alfven 
timescale ta (if ta ~ 2Ry/4:7T{p) / Bq, where (p) is the aver- 
age rest-mass density and R the stellar radius, 3 ms for 
our fiducial model); (in) the (exponential) growth rate of the 
instability should scale linearly with magnetic field strength. 
As we will discuss, all of these expectations are met. 

Figure [T] provides snapshots of the development of the in- 
stability in our fiducial star, showing projections on the (x, 
plane (upper row) and (x^y) plane (lower row) of the sim- 
ulation at times t = 1,3, 10 ms (left to right), respectively. 
These correspond to early, mid and late stages of the evolu- 
tion. Shown with vectors lines are the (global) magnetic-field 
lines, while the colors show the intensity of the toroidal mag- 
netic field only; also reported are the iso-density contours of 
the rest-mass density near the stellar surface. 

As expected, the instability develops around the neutral line 
(left column), rapidly generating a toroidal magnetic field 
in this region (middle column). The growth of this compo- 
nent continues until it reaches a comparable strength to the 
poloidal one. At this point, the growth proceeds much more 
slowly, and the magnetic field evolution is less dramatic (right 
column), as the star evolves towards a new equilibrium. As 
revealed by the different panels in Fig.[T] the development of 
the instability breaks the axisymmetry of the initial config- 
uration, leading to a complex structure with high azimuthal 
wave numbers (the m = 10 component is dominant in the 
middle column), which is eventually replaced by an m = 2 
geometry at later stages (right column). While in this complex 
evolution the (small) m = 4 component is probably inherited 
from the Cartesian coordinate system, it is interesting that the 
high-m modes develop despite the initial perturbation being 
an m = 2 one. Note also that the toroidal magnetic field pro- 
duced by the instability is concentrated in vortices (smaller at 
early times and larger at later times) and that it changes sign 
both on meridian planes (see upper row) and on the equatorial 
one (see lower row). This late-time structure is different from 
the typ ical (axisymmetric) twisted-torus discussed in previous 
works ( |Braithwaite|2QQ9l|Ciolfi et al.|2009[ [Lander & Jones] 
|2009|). We believe that this loss of axisymmetry follo ws from 
tEe conservation of magnetic helicity in ideal MHD fWoltjerl 
|T95X). Since purely poloidal (or purely toroidal) magnetic 
fields have zero helicity and the latter has to be conserved 
during a transformation in the ideal-MHD limit, the genera- 
tion of a vortex structure represents the natural way in which a 
newly-generated toroidal field will not violate the initial zero- 
helicity of the system. 

A few additional remarks are worth making about Fig. [T] 
The first one is about the evolution of the magnetic field m 
the regions right outside the star, which is essentially con- 
trolled by the resistivity there. Although the reference value 
used, tjo/Mq = 0.12, is rather high and responsible for a 
considerable decay of the magnetic field, it also allows for a 
smooth evolution and removes the development of the discon- 
tinuities which would appear in the ideal-MHD limit. As we 
will discuss later on, the qualitative behaviour of the insta- 
bility is not affected by the value of the resistivity or by the 
initial magnetic field strength. The second remark is about the 




Figure 2. Evolution of the energy in the toroidal Em, tor (red solid line) 
and poloidal £^m,poi components (black solid line), normalised to the initial 
value of the magnetic energy. The dashed lines refer to the corresponding 
evolution with 770 = 0.06 and show a smaller dissipation of the magnetic 
field (see inset) but the same dynamics in the instability. 



dynamics of the toroidal magnetic field that, as it grows to be- 
come locally comparable with the poloidal one, it also moves 
towards the stellar surface, where it can induce outflows of 
matter for smaller values of 7^0 • While a more detailed discus- 
sion of this process will be presented in a subsequent work, it 
is worth mentioning here that these winds could eject consid- 
erable amounts of matter (i.e. ~ 10~^ — 10 ~^ Mq) and thus 
have direct connections with the magnetar phenomenology. 

Figure |2] provides a more detailed analysis of the instabil- 
ity's behaviour, by monitoring the evolution of the energy in 
the toroidal £^m,tor (red solid line) and poloidal £^m,poi com- 
ponents (black solid line), normalised to the initial value of 
the magnetic energy. Initially the toroidal component is very 
small, growing slowly until a time of 2 ms. At this point, 
which is close to our Alfven-timescale estimate of 3 ms, there 
is a sudden exponential growth, which lasts for another 1 ms. 
The exponential growth, which nicely matches the predic- 
tions of the line ar-perturbation regime described by Ma rkey| 
& Tayler ( |1973 ), lasts for ~ O.Sr^ and then ceases, leaving a 
configuration which is seen to be roughly unchanged for sev- 
eral more Alfven timescales. This suggests that the configura- 
tion shown in the right column of Fig. [T] is in a quasi equilib- 
rium and is no longer susceptible to the instability. Note that 
whilst the local maxima of the field components are likely to 
dictate the system's stability, the energy eventually present in 
the toroidal-field component is only ~ 3% of the total mag- 
netic energy. This value may increase slightly on a much 
longer timescale. 

Also reported in Fig.|2]as dashed lines are the correspond- 
ing evolution of the magnetic energies when a smaller resis- 
tivity of r]o/MQ = 0.06 is used. Since the evolution of the 
instability in this case is qualitatively very similar (cf. the evo- 
lution of £^m,tor), we have confidence that our prescription 
for the resistive behaviour of the magnetic field near the stel- 
lar surface does not influence the dynamics of the instability. 
At the same time, however, a smaller resistivity is also re- 
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Figure 3. Instability growth-rate as a function of the initial magnetic field 
strength for the NSs considered (red empty circles). Note the scaling remains 
approximately linear up to — 7 x 10-*^^ G. Marked with a blue star is the 
corresponding growth-rate for 770 = 0.06. 

sponsible for a smaller decay of the poloidal magnetic field 
(see inset), which is considerably dissipated by the end of the 
simulation. While this behaviour is inevitable in a resistive 



context and has been reported also by other authors ( Braith- 
waite|2007) ), it represents an aspect of these evolutions which 



could b e improved with a fu lly consistent resistive MHD ap- 
proach (Palen zuela et al/ 2009 ). 

Another important confirmation of the perturbative analy- 
sis is offered in Fig. B] where we show the inverse of the 
growth- time r, defineathrough the exponential growth of the 
toroidal component, versus the initial magnetic-field strength 
(red empty circles). Note that the scaling is essentially linear 
for ^0 ^ 7 X 10^^ G, deviating from this for higher values, be- 
cause of the stronger magnetic tension. More specifically, the 
stronger Lorentz force will tend to oppose the fluid motions 
in the polar direction near the neutral line and which trigger 
the instability. The presence of a linear scaling is essential to 
extend our results to typical pulsar magnetic-field strengths, 
thus estimating a growth-time of ~ 10 s for a neutron star with 
Bq = 10^^ G. Also marked in Fig. [3] (blue star) is the inverse 
growth-time for the fiducial star evolved with the smaller re- 
sistivity of ijo/Mq = 0.06; again, the close similarity in the 
timescales confirms our expectation that the instability is not 
influenced by the choice of the resistivity. 

The final discussion is reserved for the potential GW signal 
emitted during the development of the instability. In Fig. |4] 
we report the GW strain in the + and x polarizations as com- 
puted from the Newtonian quadrupole formula. It is quite ap- 
parent that the signal is not of a burst type but, rather, that 
the main effect of the instability is that of triggering large- 
amplitude oscillations of the star in its fundamental F-mode. 

These GWs start emerging from the numerical noise al- 
ready at ~ 3.5 ms, but are associated to high-m oscillations 
and hence not efficient sources of GWs. However, as the 
magnetic field starts to approach the final m = 2 configu- 
ration at ~ 7 ms, the oscillations become more efficient in 
producing a GW signal (Note that 2i m = N pertubation 



Figure 4. GW strain in the + and x polarizations. Note the instability trig- 
gers large-amplitude F-mode oscillations. 

in the magnetic field leads to a m = 2N perturbation in 
the density). Because these oscillations will have a rather 
narrow spectral distribution peaked around the F-mode fre- 
quency (which is not significantly affected by the presence of 
magnetic fields), they represent very good sources of a peri- 
odic signal, potentially detectable by future advanced detec- 
tors. Defining the root- sum- square amplitude of the cross po- 

,1/2 



larization as h.^. = 



and assuming that 



the oscillations will persist undamped for :^ 0.1 — 1 s, we es- 
timate /irss = (0.54 — 1.7) X 10~^^ for a source at lOkpc. 
The corresponding signal-to-noise ratio for a detector such as 
advanced-LIGO or advanced- Virgo is S/N c:^ 1.6 — 5, thus 
potentially observable. A more detailed analysis of the spec- 
tral properties of the GW signal will be presented in a future 
work. These waveforms represent the first estimate of the con- 
version of the kinetic energy generated through the instability 
into GWs. For weaker magnetic fields, perturba tive analyses 
have suggested this coupling is much weaker ( Levin & van 
Hoven 2011 ), but more work is needed to investigate nonlin- 
early this regime. 

4. SUMMARY 

We report on numerical evolutions of the instability of 
poloidal magnetic fields in relativistic stars and the subse- 
quent generation of a mixed-field configuration in quasi- 
equilibrium. In agreement with the expectations from analytic 
perturbative studies (Markey & Tayler 1973; Wright 1973), 
we show that the instability appears after about an Alfven 
timescale, localised in the region of closed field lines. In addi- 
tion, the growth-rate of the instability has the expected linear 
scaling with magnetic field up to very large field strengths. 

The instability leads to the generation of a toroidal mag- 
netic field, which starts with a complex and high-m azimuthal 
structure, and is produced in local vortices. The instability 
ceases when the two field components are locally comparable 
in strength and the resulting configuration relaxes towards a 
geometry dominated by an m = 2 toroidal component. This 
non-axisymmetric configuration is different from the simpler 
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axisymmetric ones of earlier studies ((Ciolfi et al.|2QQ9| Lan 
|der & Jones|2009l|Braithwaite|2QQ9t[CioIti et ai.|yolO| ). How 



ever, it also represents the natural evolution of a system which 
is required to conserve its initial zero helicity in the ideal- 
MHD limit. 

The total magnetic energy of the final configuration has 
around a 3% contribution from the toroida l component ; this is 



co mparable wi th the equilibria studied by Oolfi et al.| ( |20()9] ) 

d [Lander & Jones (2009 ), but considerao^ 
of IBraithwaite (2009 ), where stratified main-sequence stars 



and [Lander & Jones ( 2009 ), but considerably lower than that 



were considered. Interestingly, the development of the in- 
stability also leads to the buoyancy of the newly-generated 
toroidal magnetic field and this can result into a mass outflow 
near the stellar surface. Additional work is needed to estab- 
lish how these winds can be related to the phenomenology 
observed in giant flares of magnetars. 

A consequence of the instability for the very strong 
magnetic fields considered here is that of triggering large- 
amplitude oscillations of the star in its F-mode. Our simu- 
lations show that the small- amplitude oscillations of the star 
are amplified by about an order of magnitude by the time the 
instability has saturated. Because these oscillations will have 
frequencies around the F-mode frequency, it is reasonable to 
perform searches for periodic signals associated to giant flares 
in magnetars. For a source at lOkpc with oscillations persist- 
ing undamped for 0.1 — 1 s, the root-sum-square amplitude 
at 1500 Hz win be /irss = (0.54 - 1.7) x 10"^^ thus leading 
to a signal-to-noise ratio S/N 1.6 — 5 for a detector such 
as advanced-LIGO or advanced- Virgo. 

During the completion of this work we have become awar e 
of a very similar analysis carried out by Lasky et al.l(|2011'), 
where the instability of purely poloidal magnetic fields in rel- 
ativistic stars was also presented. D espite the different n umer- 
ical setups, the results obtained by Lasky et al. ( [2011 ) about 
the development of the Tayler instability are in good agree- 
ment with those presented here, thus validating each other's 
conclusions. 

We are grateful to Bruno Giacomazzo, Ian Jones and Yuri 
Levin for useful comments, and Roberto De Pietri for help in 
the GW estimates. Support comes also from the "Delia Ric- 
cia" Foundation, from "CompStar", a Research Networking 
Programme of the European Science Foundation, and from 
the DFG grant SFB/Transregio 7. 
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